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available at the end of the article infected people. Various approaches have been implemented recently to study the 
dynamics of this novel infection. Mathematical models are one of the effective tools in 
this regard to understand the transmission patterns of COVID-19. In the present paper, 


we formulate a fractional epidemic model in the Caputo sense with the consideration 
of quarantine, isolation, and environmental impacts to examine the dynamics of the 
COVID-19 outbreak. The fractional models are quite useful for understanding better 
the disease epidemics as well as capture the memory and nonlocality effects. First, we 
construct the model in ordinary differential equations and further consider the 
Caputo operator to formulate its fractional derivative. We present some of the 
necessary mathematical analysis for the fractional model. Furthermore, the model is 
fitted to the reported cases in Pakistan, one of the epicenters of COVID-19 in Asia. The 
estimated value of the important threshold parameter of the model, known as the 
basic reproduction number, is evaluated theoretically and numerically. Based on the 
real fitted parameters, we obtained Ro 7 1.50. Finally, an efficient numerical scheme 
of Adams-Moulton type is used in order to simulate the fractional model. The impact 
of some of the key model parameters on the disease dynamics and its elimination are 
shown graphically for various values of noninteger order of the Caputo derivative. We 
conclude that the use of fractional epidemic model provides a better understanding 
and biologically more insights about the disease dynamics. 
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1 Introduction 

The novel virus (2019-nCoV) that is highly transmissible and pathogenic was first identi- 
fied from a single individual in Wuhan city in China. This novel infection causes a severe 
acute respiratory syndrome and it has spread across the world. The reported COVID- 
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19 confirmed cases are over 10.27 million, and there have been more than 0.5 million 
deaths till 30 June 2020 globally so far [1]. The worst affected regions due to coronavirus 
are America, Europe, Africa, South-East Asia, Western Pacific, Eastern, and Mediter- 
ranean. The initial symptoms of a COVID-19 infection include dry cough, fever, fatigue, 
and breath shortening that appear in 2-10 days and further cause pneumonia, SARS, kid- 
ney failure, and even death [2]. The pandemic has continuously spread due to absence 
of vaccine and antiviral treatments. Thus WHO announced it a global issue. The pol- 
icy makers have implemented the non-pharmaceutical intervention like social distancing, 
self-quarantine, isolation of infected, wearing mask, protective kits for medical personnel, 
and travel restrictions to minimize the disease incidence. It is also a challenging problem 
for scientists and virologist evaluating potential treatments based on ongoing clinical tri- 
als. 

Researchers suggested many mathematical models to analyze the dynamical behavior 
and spread of the novel virus which can help to predict the future situation and even 
control of the COVID-19 pandemic [3]. In the analysis of mathematical models of coro- 
navirus, the reproductive number has a significant role in describing the nonlinear dy- 
namics of physical and biological engineering problems. The reproduction number indi- 
cates that COVID-19 has been continuously increasing or has been controlled. In Pak- 
istan, 209,337 confirmed infected cases have been reported and about 4304 have lost 
their lives out of over 220 million population to date [4, 5]. The first case was reported 
in Karachi on 26 February 2020, and day by day situation is getting worse and virus is 
spreading quickly due to limiting testing. The government is unable to maintain strict 
lockdown and has imposed a smart lockdown by easing restrictions due to severe eco- 
nomic hardships, especially for labor community who earns for living to survive every 
day. To study the dynamics of COVID-19 transmission pattern, many mathematical mod- 
els provide more insight on how to control the disease spread to health authorities [6—8]. 
Fanelli and Piazza [9] studied a novel compartmental model describing the transmission 
patterns of COVID-19 in three highly infected countries. The dynamics of COVID-19 
with an impact of non-pharmaceutical interventions was studied by Ullah and Khan [10] 
on Pakistani data. The fractional mathematical models rendering the natural fact in a sys- 
tematic way as in [11, 12] and [13, 14] are used to simulate the transmission of coron- 
avirus. Different mathematical models with an effect of nonlocality and fading memory 
process by using differential operators have been presented [14, 15]. The fractional order 
epidemic models are more helpful and reliable in analyzing the dynamics of an infectious 
disease than the classical integer order models [16, 17]. The fractional order models for 
different diseases show cooperatively better fit to the real data. In [18, 19] a different frac- 
tional operator is suggested, and applications of these fractional operators are found in 
[20, 21]. Recently a Caputo fractional order COVID-19 model has been studied in [22]. 
Some other fractional mathematical models for the investigations of infectious diseases 
have been studied in [15, 23-26]. For example, the fractional diffusion equations and their 
analysis are studied in [15]. A new scheme to solve numerically the fractional order dif- 
ferential equations is utilized in [23]. New numerical investigations for the fluid in non- 
conventional media are suggested in [24]. Coronavirus modelings, simulations, and their 
possible control through a mathematical model are studied by the author in [25]. The 


spread of coronavirus in South Africa and Turkey with detailed statistical and mathemati- 
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cal results is studied in [26]. Recently, the authors have studied the analysis of coronavirus 
model in fractional derivative [27]. The use of quarantine and isolations in the model- 
ing of coronavirus is investigated in [28]. A mathematical model for the dynamical anal- 
ysis of coronavirus and its control analysis is studied by the authors in [10]. The notified 
cases of coronavirus in Saudi Arabia through a mathematical model are considered in 
[29], where the authors provide suggestions on possible controls based on the parame- 
ters. 

Environmental viral load plays an essential role in the disease incidence and is consid- 
ered to be one of the main transmission routes of COVID-19. In this study, we reformulate 
the model [28] with the impact of quarantine, isolation, and environmental effects on the 
transmission dynamics of coronavirus with the application of Caputo derivative. The pa- 
rameter values are estimated from the cumulative COVID-19 cases reported in Pakistan. 
The fractional order models provide better understanding and give more insights about 
the pandemic. The rest of the work is arranged as follows: In Sect. 2 basics preliminaries 
are presented, while the model formulation for integer case with parameter estimation and 
curve fitting is presented in Sect. 3. Model derivation and basics properties are presented 
in Sects. 4 and 5, respectively. In Sect. 6 we present the analysis of the model, while the 
numerical simulations are depicted in Sect. 7. Brief concluding remarks are presented in 
Sect. 8. 


2 Preliminaries on fractional derivative 


In order to proceed, first we recall some basic definitions regarding fractional calculus. 


Definition 1 The fractional order derivative in the Caputo case with order o for a function 
g € C" is defined as follows [19]: 


“Df (g(t)) 


t n = n-a 
1 pes s) dç, n-1<a<neN. (1) 
0 


 T(n-a) (t-s) 
Clearly, “D? (g(t)) tends to g'(t) as « > 1. 


Definition 2 The corresponding integral with order o > 0 is defined as follows: 


a 1 glot- 
I? (gt) = zl t-a] ds, O«a«Lt»0. (2) 


Definition 3 The Atangana-Baleanu-Caputo (ABC) fractional operator of order a € 
[0, 1] is defined as follows [18]: 


DeO) - POO S goel - 9" Jac, 6) 


Definition 4 The associated ABC fractional integral of order o is defined as follows: 


ABC ja _l-a d ‘ _ya-1 
Ie (eO) = Feld) + pom | (€ de, aeo w 
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Definition 5 The constant point z* is an equilibrium point of the Caputo-fractional 
model, then [30] 


CpHa(t) = g(t,z(t)), o € (0,1), (5) 
if and only if g(t, z*) = 0. 


3 The classical integer order model formulation 

In this section, we briefly discuss the integer order model of the dynamics of COVID-19 
with quarantine, isolation, and environmental load which is mainly studied for the Chi- 
nese population data in the case of fractal-fractional Atangana—Baleanu derivative [28]. 
For the mathematical model formulation, the net population N(é) at time f is further di- 
vided into mutually-exclusive sub-compartments as susceptible or healthy S(t), exposed 
E(t), infected with clinical symptoms I(t), asymptomatically infected A(t), quarantined 
class Q(t), hospitalized H(t), and the recovered R(t) individuals. The class M(t) denotes 
the environmental viral load due to the infected people. It is to be that the asymptomatic 
infected individuals are also capable of transmitting the infection. The nonlinear ordinary 
differential equations governed by these assumptions are described as follows: 


B = A — 8S(t) —AS(t), 

a - AS(t) — ((1— q)o + qp +6 + &)E(t), 
2 -(1- g)oE(t) - (c9 ci + T(t), 
45 = opE(t) — (o + 5)A(t), 


dQ (6) 
^" = e1 E(t) - (6 + m + £2)Q(t), 
22 = cI(t) + e3Q(t) - (8+ m + c3)H(t), 
R — o I(t) + A(t) + MQ) + mH() — R(E), 
eM =m I+ mA (t) = m3M(t), 
where 
y= SNE aM. (7) 


N 


subject to the following nonnegative initial conditions: 


S(0) = So > 0, E(0) = Eo > 0, I(0) = Io > 0, 
A(0) = Ao 7 0, Q(0) = Qo > 0, (8) 
H(0) 2 Ho > 0, R(0) = Ro > 0, M(0) = Mo > 0. 


The biological description of the parameters involved in COVID-19 model (6) is given in 
Table 1. 


3.1 Parameter estimation 
In this section, we estimate the model parameters with the help of the well-known sta- 
tistical technique known as nonlinear least square curve fitting approach. The confirmed 
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Table 1 Description of the model parameters and their estimated and fitted values 


Parameter Description Value (in days) Reference 
A Recruitment rate ô x N(O) Estimated 
ó atural death rate 1/(67.7 x 365) 4] 

wo ncubation period 0.1897 4] 

p ncubation period 0.1305 Fitted 

on nfected recovery rate 0.5253 Fitted 

02 Asymptomatic recovery rate 0.385 Fitted 

UE Quarantined recovery rate 0.6715 Fitted 

M Recovery rate of the hospitalized 0.6052 Fitted 

T Rate of moving from / to H class 0.306 Fitted 

£i Quarantine rate of exposed individuals 0.2395 Fitted 

E2 Hospitalization rate of Q individuals 0.310 Fitted 

e nfected disease death rate 0.0222 Fitted 

S2 Hospitalized disease death rate 0.066 Fitted 

oj Contact rate 0.6104 Fitted 

[e Disease transmission coefficient 12 x 107 Fitted 

mi Viral contribution to M by / class 0.2157 Fitted 

m2 Viral contribution to M by A class 0.2276 Fitted 

m3 Removal rate of virus from M 0.2276 Fitted 

yr Transmissibility multiple 0.5856 Fitted 


x10 
25 T T 


-9- Incidence Data 
——Model fit 


Cumulative COVID-19 Cases in Pakistan 


t-days 


Figure 1 Curve fitting (solid red line) to the confirmed infected cases using model (6), while the data are 
consider from 1 March to 30 June 2020 


infected cases from 1 March to 30 June 2020 in Pakistan are used to estimate the model 
parameter values [31]. The recruitment rate and the natural mortality rate are estimated 
from literature, and other parameters are fitted from real data. The updated reproduction 
number evaluated by using the estimated and fitted parameters is Ro ~ 1.50. The pre- 
dicted curve having a better agreement to the actual reported cases is depicted in Fig. 1. 
The model parameter values are shown in Table 1. The initial conditions used in data fit- 
ting as well as in numerical results are given as S(0) = 220,870,336, E(0) = 20,000, 7(0) = 
4, A(0) = 200, Q(0) = H(0) = M(0) = 0. 


4 Model derivation in Caputo operator 

We reformulate the COVID-19 fractional order (6) model using a Caputo fractional 
derivative in order to observe the memory effects and gain more insights about the pan- 
demic. System (6) in terms of integral form is followed by substituting the value of kernel 


as a power-law correlation function. After applying the Caputo fractional derivative of 
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order a — 1, we get 


Cpe-t[65] = €pa-tz. © YA — 850) — AS), 

Cpe 42] = CDO DASE) - (0 — yo + ep + ô + EEO), 

Cpe-i[d] = Cpe-1r (1 — g)e(t) - (e + ô+ «i + TOL, 

cpe-144] = CDr- [pp E(t) — (02 + 5)A(0)], és 
cpe-1[49] = Cpe le Elt) - (8 + m + e) Q(0], 

Cpe-1[84] = Cpe Miri (t)  e3Q(t) — (8 + m  c2)H(t], 

Cpe B] = Dey lolt) + Alt) + m Q(t) + mH(t) - RE), 

Cpe-1[4M] = Cpe [mI + mA) - m3M(t)]. 


Both are the inverse operators, we get the COVID-19 model in the Caputo operator form 


CD*'S(t) = A — S(t) - AS(t), 

CDYE(t) = AS(t) - (0 — q)o + qp + 8 + &)E(t), 
CD?I(t) = (1 — e)oE(t) - (o1 + 8 + c + T(E), 
CD? A(t) = gpE(t) — (o» + 8)A (t), 


(10) 

CD? Q(t) = eE(t) - (8 + m + £2)Q(t), 

CD?H(t) = TI(t) + eQ() — (8 + m + c3)H(, 

CDI R(t) = I(t) + o2A(t) + m Q(t) + nH (t) - 5R(t), 

CD? M(t) = mi + m2A(t) —m3M(t). 

The initial conditions are as follows: 
S(0) = 0, E(0) = 0, I(0) = 0, A(0) 7 0, ii 

11 


Q(0) = 0, H(0)>0, R(0) > 0, M(0) = 0. 


5 Basic properties 

5.1 Invariant region 

The dynamics of Caputo fractional model (10) is explored in a feasible region Q C IR? such 
that Q = ((S(£), E(t), I(t), A(t), Q(t), H(t), R(t) € R7 : N(t) x $,M(t) € R; :< 5an), 
Lemma 1 The region Q C R? is positively invariant with nonnegative initial conditions 
for model (10) in R7. 


Proof After summing the components of human population in model (10), we obtain a 
total population as follows: 


CD*N(t) = CD S(t) + CD E(t) + CD I(t) + -D* A(t) + CD* Q(t) 


+ €D*H(t) + -D?R(t), 
and then we have 


CD*N(t) + àN(t) < A. 
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After applying the Laplace transform, we have 


a-1 


A 
Ns) < s(s* + à) ir sS" 


By considering the inverse Laplace, we arrive at 
N(t) < N(0)Eq,1 (8t*) + At” E, a1 (6t"), 


where the Mittag-Leffler function is shown by 


eo n 


y 
E = y c —X Ó 
pG) = I (on + B) 


and the Laplace transform is 


s£ 


LP Engla) = 2—. 


Thus, the solution of the model with the nonnegative conditions in Q remains in Q. So, 


the region Q is positively invariant and attracts all the solutions in R^. 


Now, for the positivity of the system solution, let 
7 7 T 
Ri ={yeR’|y=0} and y(t) = (SW), EQ).I(),A(0, QU), HORO) . 


Corollary 1 [32] Suppose that g(t) € C[m, n] and ©D%g(t) € (m,n],where a € (0,1). Then 
if 

(i). CDe'e(t) >0, Vye(m,n)then g(t) is nondecreasing; 

(ii). CD? g(t) <0, Vy e (m,n), then g(t) is nonincreasing. 
5.2 Positivity and boundedness 


Proposition 1 The solution of model (10) is nonnegative and bounded for all (S(O), E(0), 
1(0), A(0), Q(0), H(0), R(0)) € R7 for t > 0. 


Proof To show that the solution of the model is nonnegative, it is required to show that on 
each hyper-plane bounding the positive orthant there is a the vector field point R7. From 
system (10), we have 


CD*S|so- A»0, —|CD*E|go- A$- 0,  CD*I |r-0= (1— g)oE > 0, 
“DA |4-0= GPE>0, | €*D?Ql|g-o-&E-0, ©DYH |u-o- TI + €3Q > 0, 


CD*R |g-o7 O11 + 03A + MQ + noH > 0. 
The solution is 


N(t) < N(0)E, 1 (-0t") + At” Ey. (t^). 
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Also, the Mittag-Leffler function is bounded V t > 0. Therefore, we have 


e| - 


lim N(t) < 


The solution of the system will remain in R7 by using Corollary 1, and hence a biologically 


feasible region is constructed as follows: 


Q = {(S,E,LA,Q,H,R) € RÌ : S,E,I,A,Q,H,R = 0}. 


Since all the terms are positive, the solution of model (10) is bounded. 


6 Analysis of the model 
6.1 Disease-free equilibrium (DFE) 


For the equilibrium points of fractional model (10), we have 


©D* S(t) = €D*E(t) = CD(t) = CD% A(t) = DF Q(t) = £D? H(t) = ČD R(t) 


= “DY M(t) = 0. 


The disease-free equilibrium is denoted by E; = (S9, E9, I9, A9, Q9, H9, R°, M?) and is given 
by 


A 
EO = (50000000). (12) 


6.2 The reproduction number 
The next generation technique is used to derive the reproduction number Ro for the dy- 
namics of disease. Let x = (E, I, A, Q, H, M)T, then we have 


dx 
— -E-V, 13 
3i (13) 
where 
0 va 0 22 0 
0 0 0 0 0 0 
jo 0000 
loo o o o Of’ 
0 0 0 0 0 0 
0 0 0 0 0 0 
Ko 0 0 0 0 0 
(p-1)o k 0 0 0 0 
v= —Qp 0 k 0 0 0 
—Ói 0 0 k 0 0 
0 -T 0 —& k4 O0 
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The next generation matrix is of the form 


y kokakamao(1—9) A y NG AC 
SED : EDU TS +d m * Tom) ^u * bi 0 0 "o 
0 0 0 0 0 0 
rv = 0 0 0 O 
0 0 0 0 0 0 
0 0 0 0 0 0 
0 0 0 0 0 0 
where 
ko = (1- g)o + gp + ô+ £1, ky =o, +6+6,+T, ky = 07 + ô, 
k3=5 +) t£», ky =8 + 2+ 6», 


(Kykakam»qp + kokakamie — kokskamyeog) ^ çz 


d 
Koki kok3kqgm36 


The spectral radius of the next generation matrix is shown by 


_ kidmshigpy + kxdm30b1 - kotypo Am, 


R 
" kokikom3ó 


E Kitopo Nm» + katyo Nm, — koótomao 


kok kom36 


6.3 Local stability of DFE 


(14) 


Theorem 1 For any two positive integers q and r, let gcd(q,r) = 1 for a = 4 and M = n, the 


disease-free equilibrium is locally asymptotically stable if | arg(A)| > 55; for all roots of the 


associated characteristic equation 


deitding ee E] —J(E})) = 0. 


(15) 


Proof For local stability, the linearized system is the Jacobian of system at disease-free 


state as follows: 


- 0 -ü -va 0 0 0 Œ 

0 -h à wh 0 0 0 % 

0 (1-g)o -k 0 0 0 0 0 

(ES) = 0 9) o -k 0 0 0 0 
s 0 £ 0 0 -k 0 0 0 

0 0 T 0 à -k 0 0 

(0) (0) 01 02 nı n —ô (0) 

(0) (0) mı ma (0) (0) 0 -ma 


The evaluation of determinant equation (15) implies 


(At + 8) Q4 + k3)(A4 + k4) (A? + a A*4 + aa?i a3 + 44) - 0, 


(16) 


(17) 
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the eigenvalues —ô, —ô, —k3, —k4 have negative real parts and the others can be found from 
the last factor of (17). The co-efficient is of the form 
a, = ko + kj + k + mg, 
az = koka(1 — R4) + (ko + ky + ko)m3 + kiko, 
43 = kokika(1 — R2) + kokoma(1 — R3) + kKokyma(1 — Ra) + kıikomz 
- & (ekipy * m(ppY + (1— 9)o)), 
a, = kokykoma(1 — Ro). 


Clearly, a; for i = 1,2,...,4 are all positive if Ro < 1. The arguments of the roots of equa- 
tions (A4 + 8)? = 0, (A4 + kg) = O and (A + ka) = O are similar, that is, 


2m m T 
> > A 
m M 2 


arg) = » +k where k = 0,1,2,..., (m — 1). 


In a similar fashion, we also find that the arguments of the equation (A4 + 21474 + 45474 + 
43A t a4) = O are all greater than aM if Ro < 1, having an argument less than aM for Ro > 1. 


The DFE is locally asymptotically stable for Ro « 1. 


6.4 Global stability of DFE 
The Lyapunov function approach is used to proceed to the result for the GAS of the pro- 
posed model at the disease-free state. For this, we have the following theorem. 


Theorem 2 The DFE of the Caputo COVID-19 model is GAS if Ro « 1. 
Proof Consider the following appropriate Lyapunov function: 
F) = YE + Yor + Y3A + YAM, (18) 


where the coefficients Yj, for j = 1, 2,3, 4, are unknown positive constants, and they will be 
chosen later. The Caputo-fractional derivative of $(£), along model (10), yields 
CD*$(t) = Yı DE + Yo°DYI + Y3°D%A + Y4“ D?M, (19) 
“DY F(t) = Y1 (AS(t) — KoE(t)) + Yo((1 — g)@E(t) — ki I(t)) 
+ Y3(gpE(t) — koA(t)) + Ya(mi(t) + mA(t) — m3M(t)) 
< (Yit - Yok, + Yam (0) + (Yow — Yk + Yam3)A(t) 
+ (Y2(1— 9)o — Yiko + Ysgp)E(t) + (£58? — Yaz) M(t), S < N, S 
€ (X46 — Yakı + Yam I(t) + (Vidi — Yako + Yam3)A(t) 


Y3(1 - g)w + Yaop 
Yı ko 


+ zi 1)so + (S? - Yam3)M(t). 


Let us choose 


óCq- A 
Yı = maó, Y, = ee Ys = m3yl16 + Amo, Y, = A6». 


k 1 
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CD* (t) < kom35(Ro — WE. (20) 


Hence it follows that CD? < 0 for Ro < 1, all the parameters and variables are nonnega- 
tive with CD" = 0 iff E = I = A = M = 0. Thus (E,L, A, M) — (0,0,0,0) as t > oo. By using 
E=I=A =M =O inthe first and second last three equations of (10) implies that S > 4, 


and Q, H, R > 0ast — oo. Thus, using Lyapunov stability theorems for the fractional case 


developed in [30], the solution of Caputo model (10) with nonnegative initial conditions 
approaches to Ej as t — oo ina feasible region. Thus, it follows that the DFE of model (10) 
is GAS. 


6.5 Existence of endemic equilibrium point 
We denote the endemic equilibrium for fractional coronavirus model given in (10) by Eoy = 
(S'" EM p AS Q"*, H**, R**, M**) such that 


St = sts 
3k OK 
E** = À m ; 
pe = OOE _ -pors 
~ kı ~ koky , 
A** = gpE** _ gpa**s** 
T k ^ kok ? (21) 
Q** _ gET _ e*s 
BE GT 
H* = tI**+e2Q** _ rKa(1-9)omA** S** ck1e1e3A** S** 
k4 kokikok3kq J 
M** m2A** +m I** kym»gpX** S** «komi(1-g)oX** S** 
i ma ~ kok ko ma : 
R* oA «ng H** +011** «gi Q** 
- 5 . 
Inserting the above result into the force of infection state, we get 
T* 4 WA 
ans d E E dM. (22) 


N* 


Substituting (21) into (22) shows that the nonzero equilibria of the model satisfy the 


quadratic equation 
bod**? + b1A** + by = 0, (23) 
where 


bo = kokikams (ko (eia (e2(8 + n2) + k4 (ô + m)) + ksks) + gokıkzk4p(ô + &2)), 
by = kokık28mz (ko (e1kı (e2(6 + n2) + Ka(8 + 1)) + Kakz) + økıkzkap (8 — Wei + 02)) 
+ ok A((1 — p)kgmo — vk, oma) + KK skim, 


b; = Kika Ka, ms(1— Ro). 
Also, we have 


ks = 1(5 + no(1 — g)o + ka((1 — 9)o(8 + m) + a6), 


ks = ko(- ko (ela (e(8 + m2) + Ka (8 1) — Ksks) — elakskap(8 + 02), 
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k7 = t(0 + m(1 — g)o + ka((1 — g)o(8 + m — &) + ki). 


Here, it should be noted that the coefficient of (23), bo > 0, and b» is positive (negative) 
Ro < 1, (>). The following result is established. 


Theorem 3 The proposed model has: 
(i) a unique (Ej1) endemic equilibrium if b; < 0 & Ro > 1; 
(ii) unique E** if by « 0, ^b; = 0) v b? — 4bob = 0; 
(ii) if bi « 0, bz > 0 and the discriminant is greater than zero, then there exist two 


endemic equilibria. 


(iv) Otherwise, no endemic equilibrium exists. 


Thus from case (i) of Theorem 3, it is clear that model (10) has a unique positive endemic 
equilibrium when 7, > 1, and point (iii) indicates the possibility of backward bifurcation 
when Fo « 1. 


7 Numerical results 

This section presents simulation and discussion for the Caputo COVID-19 model (10). 
The proposed fractional model is solved numerically using a generalized predictor- 
corrector of the Adams-Bashforth-Moulton method [33, 34]. The biological parame- 
ters estimated from the actual data reported in Pakistan from 1 March to 30 June 2020 
and tabulated in Table 1 are utilized to obtain the simulation results. We varied different 
model parameters and the order o of the Caputo operator in order to explain the role of 
various parameters and memory index on the disease transmission patterns and control. 
The impact of contact rates ¢1, £5, V (a measure of social distancing effects), quarantine 
or contact-tracing policy £1, and hospitalization or self-isolation rate t is depicted graph- 
ically in Figs. 2-6. The impact of environmental viral load due to the symptomatic (m1) 
and asymptomatic infected individuals (775) is depicted in Figs. 7 and 8, respectively. Fi- 
nally, the role of removal rate of virus from the environment (via disinfection spray etc.) 


is shown in Fig. 9. 


7.1 Effect of reduction in contact rates (a measure of social distancing) 

First of all, we simulate COVID-19 model (10) for different values of effective contacts £1 
in order to evaluate the impact of social distancing on the dynamics of cumulative symp- 
tomatic, asymptomatic, and hospitalized infected individuals. The simulation results for 
different four values of fractional order are shown in Fig. 2(a)- (d). It can be observed from 
these interpretations that the reduction in the effective contact rates 7; up to 50 percent 
to its baseline values reduces the peak of cumulative COVID-19 infected cases very well. 
The behavior is the same, seen for all four values of fractional order, although slightly 
faster decrease in pandemic peak is observed for smaller values of o. The impact of disease 
transmission rate due to the environmental viral load ¢2 is demonstrated in Fig. 3(a)- (d). 
With the reduction in this rate, a rapid decrease is observed in the total of symptomatic, 
asymptomatic, and hospitalized infected individuals and even becomes more biologically 
significant for smaller values of fractional order. This interpretation shows that the envi- 


ronmental transmission plays a significant role in disease prevalence. Finally, the influence 


Aba Oud et al. Advances in Difference Equations 


(2021) 2021:106 


P =1 s =0.95 
Nt. o uat l , g : : 
nsns — (ED Baseline 
er "F | — c =0.4883 ] 
sl: sa | — (70-3608 | 
€ € | |—¢-024 
I, I, 
HE $; 
4r 4 
2f L4 
0 1 0 
0 50 100 150 200 250 300 0 50 100 150 200 250 300 
Time(days) Time(days) 
(a) (b) 
; =0.90 ; =0.85 
44 x10 i Q T : 4 x10 o T I 
—— (706104, Baseline — 706104, Baseline 
12 —(,=0.4883 4 12H — (,-0.4883 4 
—(,70. 3608 —ü23. 3608 
0 1 10H LC 202441 1 
2 70.2441 2 1 
r 8 r ak 
f; $ et 
4 ab 
£ : jen 
" 1 ; 1 P 7 , 1 
0 50 100 150 200 250 300 0 50 100 150 200 250 300 
Time(days) Time(days) 
(c) (d) 
Figure 2 The impact of contact rate £1 on cumulative symptomatic, asymptomatic, and hospitalized 
COVID-19 individuals for œ = 1,œ = 0.95,œ = 0.90, œ = 0.85 


"i 
"ni 


az1 


I()+A(t)+H(t) 


I(t) +A(t)+H(t) 


x10 


(c) 


0 50 100 150 200 250 300 D] 50 100 150 200 250 300 
Time(days) Time(days) 
(a) (b) 
; a=0.90 s a=0.85 
ase r T I a5 T T I 
15 «108 
ar N al 3 
(712 «10 
NI —g1 x 10° p 3 
I —— (08x 10* I 
i 2 tc 08% 4 i 2 
€ € 
b 15 1 15 
1 1 
05 05 
0 0 
0 50 100 150 200 250 300 0 50 100 150 200 250 300 
Time(days) Time(days) 


(d) 


Figure 3 The impact of contact rate € on cumulative symptomatic, asymptomatic, and hospitalized 
COVID-19 individuals for œ = 1,a@ = 0.95,œ =0.90,a@ = 0.85 


Page 13 of 19 


Aba Oud et al. Advances in Difference Equations 


(2021) 


2021:106 


ug aet T T PES T iid T T 
tap H 
tor H 
f 11 
=I | = 
Al | 
à | 
0 50 100 150 200 250 300 0 50 100 150 200 250 300 
Time(days) Time(days) 
(a) (b) 
"t o-0.90 Nt o-0.85 
— (= 0.5856, Baseline —1= 0.5856, Baseline 
12 
10 
i i 
= i =. 
2 
" i ; , 1 
0 50 100 150 200 250 300 0 50 100 150 200 250 300 
Time(days) Time(days) 
(c) (d) 
Figure 4 The impact of y on cumulative symptomatic, A, and H COVID-19 individuals for 
a=1,a =0.95,a = 0.90,a@ = 0.85 


0 EJ 100 150 


Time(days) 


(c) 


a=1,a =0.95,a = 0.90,a@ = 0.85 


200 


250 300 


a= y 00.95 
PEU : : ? "P ; ; : ' 
—,= 03101, Baseline —«,= 0.3101, Baseline 
£ £ 
T I 
i E 
€ € 
Ed E 
i d 
€ € 
0 50 100 150 200 250 300 0 50 100 150 200 250 300 
Time(days) Time(days) 
(a) (b) 
vet? a=0.90 ; az0.85 
om i T T r r ug 1 r T 
—,= 0.3101, Baseline —<,= 0.3101, Baseline! 
Ji 03721 ] 
"i z04841 | 
€ | | — «04961 2 
I I 
FAI ils 
€. | ¢ 
E t 
m J 
ab y 


Time(days) 


(d) 
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Figure 9 The impact of removal rate m3 on symptomatic, A, and H COVID-19 individuals for 
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of the reduction in y is depicted in Fig. 4(a)—(d) for four different values of a. A reason- 
able decrease in the pandemic peak is observed with a reduction of contact rate y. Over- 
all, from these interpretations, we conclude that maintaining a strict social distancing and 


avoiding public gathering are needed to control the transmission of infection in the future. 


7.2 Effect of contact-tracing and hospitalization/self-isolation 

The quarantine contact-tracing policy of exposed individuals and hospitalization/isolation 
of confirmed infected cases is another effective strategy to control the COVID-19 inci- 
dence. Here, we simulate the model with different levels of quarantine and hospitalization 
rates. The resulting graphical interpretation is shown in Figs. 5(a)- (d) and 6(a)- (d), re- 
spectively. A dramatic decrease in pandemic peak is observed with an increase (up to 50 
percent) in the contact-tracing rate as shown in Fig. 5(a)- (d). The same behavior is ob- 
tained for all values of o. The influence of the increase in hospitalization of confirmed 
infected individuals is analyzed in Fig. 6(a)—(d). These results suggest that in the absence 
of treatment or vaccine, a proper and effective contact-tracing and self-isolation policy 
should be adopted until the elimination of pandemic. 


7.3 Effect of environmental viral load 

The environmental viral load is one of the main routes of the COVID-19 infection trans- 
mission. In this section, we simulate model (10) with the reduction in the viral contribution 
to the environment by symptomatic and asymptomatic COVID-19 infective individuals at 
different rates. Figure 7(a)-(d) demonstrates the impact of the reduction in mı with dif- 
ferent levels. It is observed that a reduction in mı (up to 50 percent) reduces the infected 
curve with a reasonable rate as seen in Fig. 7(a)- (d). The decrease in pandemic peaks is 
relatively more significant for smaller values of the fractional order of Caputo operator. 
The role of m (viral contribution of asymptomatic individuals to the environment) is an- 
alyzed in Fig. 8(a)- (d). It can be seen that a reduction in m3 at 50 percent to its baseline 
values given in Table 1 reduces the total number of symptomatic, asymptomatic, and hos- 
pitalized infected individuals very well. The similar behavior is obtained for all values of 
fractional order o. Finally, the influence ofthe removal rate of virus from the environment 
is depicted in Fig. 9(a)- (d). It is observed that enhancement in mg at 40 percent to its 
baseline value dramatically reduces the infection curves as shown in Fig. 9(a)- (d). This 
interpretation suggests that proper disinfection spray is necessary to remove the viruses 


from the environment in order to reduce the disease prevalence. 


8 Conclusion 

The environmental viral load plays an essential role in the spread of the current COVID- 
19 pandemic. Limited research has been done on the environmental impacts on disease 
dynamics. In this paper, we studied the transmission dynamics of COVID-19 pandemic 
with quarantine, hospitalization (or self-isolation), and the environmental viral load im- 
pacts through a Caputo fractional model. Initially, the model is presented with classical 
integer order differential equations and then reformulated using fractional order opera- 
tor with the power-law kernel. The basic mathematical analysis of the Caputo COVID-19 
model is performed and the local and global stability of the disease-free case is proven. 
Utilizing the nonlinear least square approach, the model parameters are estimated and 
fitted for the COVID-19 actual cases recorded in Pakistan. The most important threshold 
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quantity, known as the basic reproduction number, is presented both theoretically and 
numerically. The predictor-corrector iterative scheme of Adams—Moulton type is applied 
in order to solve the model numerically. The real estimated and fitted parameters are used 
to simulate the Caputo COVID-19 model. We depicted and discussed the simulation re- 
sults for three different sets of model parameters showing the impact on the dynamics 
cumulative symptomatic, asymptomatic, and hospitalized infected cases. Firstly, we ana- 
lyzed the influence of reduction in contact rates (a measure of social distancing) on the 
infected curves. The graphical interpretation of this case revealed that with a reduction 
of effective contacts with infected individuals, a significant decrease in the cumulative in- 
fected cases is seen. It is further observed that the reduction in the environmental viral 
load dramatically reduces the pandemic peak. Secondly, we explored the impact of the 
most commonly used interventions i.e., quarantine (or contact-tracing policy) and hos- 
pitalization (or self-isolation) on cumulative infected population. These graphical inter- 
pretations demonstrated that contact-tracing of exposed people is more effective in the 
elimination of infection. It can be seen that with an increase in the quarantine rate up to 
50 percent to the base line value reduces the infected cases significantly. Finally, we have 
shown the impact of variation in the environmental viral load due to symptomatic and 
asymptomatic COVID-19 infected individuals. It is found that reducing the rate of viral 
release into the environment by asymptomatic infected individuals decreases the disease 
burden very well. It is because the asymptomatic infected individuals are unaware of the 
infection and freely continue the daily routine without any necessary safety precautions 
like the use of face mask and sanitizer. Moreover, we have also depicted the impact of 
variation in the removal rate of virus from the environment (or surfaces) on the disease 
prevalence. A dramatic reduction in the infected cases is found with an enhancement in 
this rate, which demonstrates that disinfection spray is very helpful in the disease elim- 
ination. We believe that the investigations in the paper will be beneficial for the health 
and decision-making authorities to combat the disease in the community. In the future, 
the present model can be reformulated by incorporating some suitable time-dependent 


control interventions using optimal control theory. 
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